/*
 * @(#)Simplex.java        1.0 2000/06/26
 *
 * This file is part of Xfuzzy 3.0, a design environment for fuzzy logic
 * based systems.
 *
 * (c) 2000 IMSE-CNM. The authors may be contacted by the email address:
 *                    xfuzzy-team@imse.cnm.es
 *
 * Xfuzzy is free software; you can redistribute it and/or modify it
 * under the terms of the GNU General Public License as published by
 * the Free Software Foundation.
 *
 * Xfuzzy is distributed in the hope that it will be useful, but WITHOUT
 * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
 * FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
 * for more details.
 */

package xfuzzy.xfsl.model.algorithm;

import xfuzzy.xfsl.model.*;
import xfuzzy.lang.*;
import java.util.Random;
import xfuzzy.xfds.XfdsDataSet;

/**
 * Algoritmo de descenso por Simplex
 * 
 * @author Francisco Jos� Moreno Velo
 *
 * @see "Press, W.H., Teukolsky, S.A., Vetterling, W.T., Flannery, B.P., 
 * Numerical Recipes in C, Cambridge University Press, 1997."
 */
public class Simplex extends XfslAlgorithm {

	//----------------------------------------------------------------------------//
	//                            MIEMBROS PRIVADOS                               //
	//----------------------------------------------------------------------------//

	/**
	 * Constante para generar n�imeros aleatorios
	 */
	private double pert;
	
	/**
	 * Factor de expansi�n
	 */
	private double inc;
	
	/**
	 * Factor de contracci�n
	 */
	private double dec;
	
	/**
	 * Descripci�n de las N+1 configuraciones del simplex
	 */
	private double[][] simplex;
	
	/**
	 * Suma de los valores de cada par�metro en el simplex,
	 * almacenado para simplificar los c�lculos
	 */
	private double[] simsum;
	
	/**
	 * Evaluaciones de cada configuraci�n del Simplex
	 */
	private XfslEvaluation[] simeval;
	
	/**
	 * �ndice de la mejor configuraci�n del simplex
	 */
	private int ilo;
	
	/**
	 * �ndice de la peor configuraci�n del simplex
	 */
	private int ihi;
	
	/**
	 * �ndice de la segunda peor configuraci�n del simplex
	 */
	private int inhi;
	
	/**
	 * Mejor error de la �ltima iteraci�n
	 */
	private double lasterror;
	
	/**
	 * Generador de n�meros aleatorios
	 */
	private Random random;

	//----------------------------------------------------------------------------//
	//                                CONSTRUCTOR                                 //
	//----------------------------------------------------------------------------//

	/**
	 * Constructor
	 */
	public Simplex() {
		this.pert = 1;
		this.inc = 2;
		this.dec = 0.5;
		this.random = new Random();
	}

	//----------------------------------------------------------------------------//
	//                             M�TODOS P�BLICOS                               //
	//----------------------------------------------------------------------------//

	/**
	 * Devuelve el c�digo de identificaci�n del algoritmo
	 */
	public int getCode() {
		return SIMPLEX;
	}

	/**
	 * Actualiza los par�metros de configuraci�n del algoritmo
	 */
	public void setParameters(double[] param) throws XflException {
		if(param.length != 3) throw new XflException(26);
		pert = test(param[0], POSITIVE);
		inc = test(param[1], INCREASE);
		dec = test(param[2], DECREASE);
	}

	/**
	 * Obtiene los par�metros de configuraci�n del algoritmo
	 */
	public XfslAlgorithmParam[] getParams() {
		XfslAlgorithmParam[] pp = new XfslAlgorithmParam[3];
		pp[0] = new XfslAlgorithmParam(pert, 1, POSITIVE, "Perturbation");
		pp[1] = new XfslAlgorithmParam(inc, 2, INCREASE, "Expansion Factor");
		pp[2] = new XfslAlgorithmParam(dec, 0.5, DECREASE, "Contraction Factor");
		return pp;
	}

	/**
	 * Ejecuta una iteracion del algoritmo
	 */
	public XfslEvaluation iteration(Specification spec, XfdsDataSet pattern,
			XfslErrorFunction ef) throws XflException {
		Parameter[] param = spec.getAdjustable();
		if(init) {
			init = false;
			simplex = new double[param.length+1][param.length];
			simsum = new double[param.length];
			simeval = new XfslEvaluation[param.length+1];
			for(int j=0; j<param.length; j++) simplex[0][j] = param[j].value;
			for(int j=0; j<param.length; j++) simsum[j] += param[j].value;
			simeval[0] = ef.evaluate(spec,pattern,1.0);
			lasterror = simeval[0].error;
			for(int i=1; i<simplex.length; i++) {
				perturb(spec);
				for(int j=0; j<param.length; j++) simplex[i][j] = param[j].value;
				for(int j=0; j<param.length; j++) simsum[j] += param[j].value;
				simeval[i] = ef.evaluate(spec,pattern,lasterror);
				for(int j=0; j<param.length; j++) param[j].value = simplex[0][j];
			}
			getLimits();
			for(int j=0; j<param.length; j++) param[j].value = simplex[ilo][j];
			lasterror = simeval[ilo].error;
			return simeval[ilo];
		}

		XfslEvaluation simtry = amotry(spec,pattern,ef,param,ihi,-1.0);
		if(simtry.error <= simeval[ilo].error)
			simtry = amotry(spec,pattern,ef,param,ihi,inc);
		else if(simtry.error >= simeval[inhi].error) {
			XfslEvaluation simsave = simeval[ihi];
			simtry = amotry(spec,pattern,ef,param,ihi,dec);
			if(simtry.error >= simsave.error) {
				for(int i=0; i<simplex.length; i++) if(i != ilo) {
					for(int j=0; j<param.length; j++)
						param[j].setDesp(dec*(simplex[i][j] - simplex[ilo][j]));
					spec.update();
					simeval[i] = ef.evaluate(spec,pattern,1.0);
					for(int j=0; j<param.length; j++) simplex[i][j] = param[j].value;
					for(int j=0; j<param.length; j++) param[j].value = simplex[ilo][j];
				}
				for(int j=0; j<param.length; j++) {
					simsum[j] = 0;
					for(int i=0; i<simplex.length; i++) simsum[j] += simplex[i][j];
				}
			}
		}

		getLimits();
		for(int j=0; j<param.length; j++) param[j].value = simplex[ilo][j];
		simeval[ilo].var = (lasterror - simeval[ilo].error)/lasterror;
		lasterror = simeval[ilo].error;
		return simeval[ilo];
	}

	//----------------------------------------------------------------------------//
	//                             M�TODOS PRIVADOS                               //
	//----------------------------------------------------------------------------//

	/**
	 * Evalua una proyecci�n de factor fac
	 */
	private XfslEvaluation amotry(Specification spec, XfdsDataSet pattern,
			XfslErrorFunction ef, Parameter[] param, int ihi, double fac) {
		double[] ptry = new double[ simplex[ihi].length ];
		double fac1 = (1.0 - fac)/simplex[ihi].length;
		double fac2 = fac1 - fac;
		for(int i=0; i<param.length; i++) {
			param[i].value = simplex[ihi][i];
			ptry[i] = simsum[i]*fac1 - simplex[ihi][i]*fac2; 
			param[i].setDesp(ptry[i] - simplex[ihi][i]);
		}
		spec.update();
		XfslEvaluation ytry = ef.evaluate(spec,pattern,1.0);
		if(ytry.error < simeval[ihi].error) {
			simeval[ihi] = ytry;
			for(int i=0; i<param.length; i++) {
				simsum[i] += param[i].value - simplex[ihi][i];
				simplex[ihi][i] = param[i].value;
			}
		}
		for(int i=0; i<param.length; i++) param[i].value = simplex[ilo][i];
		return ytry;
	}

	/**
	 * Actualiza los indices del mejor, peor y antepeor
	 */
	private void getLimits() {
		ilo = 0;
		ihi = (simeval[0].error > simeval[1].error ? 0 : 1);
		inhi = (simeval[0].error > simeval[1].error ? 1 : 0);
		for(int i=0; i<simplex.length; i++) {
			if(simeval[i].error <= simeval[ilo].error) ilo = i;
			if(simeval[i].error > simeval[ihi].error) { inhi = ihi; ihi = i; }
			else if(simeval[i].error > simeval[inhi].error && i!=ihi) inhi = i;
		}
	}

	/**
	 * Realiza una perturbaci�n del sistema	
	 */
	private void perturb(Specification spec) {
		Type[] type = spec.getTypes();
		for(int i=0; i<type.length; i++) {
			if(!type[i].isAdjustable()) continue;

			Family fam[] = type[i].getFamilies();
			for(int j=0; j<fam.length; j++) {
				if(!fam[j].isAdjustable()) continue;
				Parameter[] pp = fam[j].getParameters();
				double[] val = fam[j].get();
				do {
					for(int k=0; k<pp.length; k++) {
						if(!pp[k].isAdjustable()) continue;
						pp[k].value = val[k] + pert*(random.nextDouble()-0.5);
					}
				} while(!fam[j].test());
			}

			ParamMemFunc[] mf = type[i].getParamMembershipFunctions();
			for(int j=0; j<mf.length; j++) {
				if(!mf[j].isAdjustable()) continue;
				Parameter[] pp = mf[j].getParameters();
				double[] val = mf[j].get();
				do {
					for(int k=0; k<pp.length; k++) {
						if(!pp[k].isAdjustable()) continue;
						pp[k].value = val[k] + pert*(random.nextDouble()-0.5);
					}
				} while(!mf[j].test());
			}
		}
	}
}

